Exact solution of a model of qubit decoherence due to telegraph noise 
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We present a general and exact formalism for finding the evolution of a quantum system subject 
to external telegraph noise. The various qubit decoherence rates are determined by the eigenvalues 
of a transfer matrix. The formalism can be applied to a qubit subject to an arbitrary combination 
of dephasing and relaxational telegraph noise, in contrast to existing non-perturbative methods that 
treat only one or the other of these limits. As the applications: 1) We obtain the full qubit dynamics 
on time scales short compared with the enviromenal correlation times. In the strong coupling cases 
this reveals unexpected oscillations and induced magnetization components; 2) We find in strong 
coupling case strong violations of the widely used relation I/T2 = l/2Ti + l/T^, which is a result 
of perturbation theory; 3) We discuss the effects of bang-bang and spin-echo controls of the qubit 
dynamics in general settings of the telegraph noises. Finally, we discuss the extension of the method 
to the cases of many telegraph noise sources and multiple qubits. The method still works when 
white noise is also present. 



I. INTRODUCTION 

Decoherence of quantum systems is a fundamental is- 
sue with impHcations across all branches of physics. In 
this context, two-level quantum systems have served as a 
very useful paradigm. The subject of two-level systems 
in random time-dependent fields originated decades ago 
in the context of spin resonance [1]. Later experiments 
on macroscopic quantum coherence also focused atten- 
tion on this problem [2]. The recent surge of interest 
in quantum computing and quantum control has reju- 
venated the field, and much work has gone into solving 
various models of decoherence. All physical realizations 
of qubits are subject to external noise since all couple, 
however weakly, to the external environment. The most 
popular models are those that explicitly involve a bath 
whose degrees of freedom must be traced out [3]. In many 
cases, however, one can neglect the flow of quantum in- 
formation from the bath to the system, and then it is 
sufficient to consider the qubit as being subject to ran- 
dom classical external fields. The conditions under which 
this assumption is valid have been considered in detail by 
Galperin et al. [4]. Here we merely note that when the 
coupling of the bath to its thermalizing external envi- 
ronment is very strong or on time scales longer than the 
characteristic microscopic times of the bath, we expect 
that even fully quantum system-bath models reduce to 
this case. 

The most important type of noise in solid-state sys- 
tems is telegraph noise: the qubit interacts with one or 
more random fluctuators in its neighborhood, and these 
fluctuators go back and forth between only two states 
[5]. Usually the qubit interacts with many such fluctu- 
ators. Depending on the distribution of fluctuator tran- 
sition rates this may give rise to 1// or other types of 
noise. In this paper, however, we shall focus on the case 
of a single fluctuator, and will only mention the general- 
ization to many fluctuators at the end. This special case 
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FIG. 1: One realization of telegraph noise drawn from a dis- 
tribution with p = 0.1, 5 = and At — 0.1 (in arbitrary 
units). See the text for details. 

is also extremely important in experimental practice [6]. 
Thus we investigate the following model: in a sequence 
of N time intervals each of length At, the fluctuator pro- 
ducing the noise undergoes a sequence of states labelled 
by {si, S2, Sat}, Si = ±1 as shown in Fig. 1, and so 
the hamiltonian within the i-th interval is given by 

H,^ = -^{Bo + s.g) ■ a, (1) 

where <yx.y,z are the Pauli matrices acting on the qubit 
spin, Bq = Bqz pointing in the z-direction for definite- 
ness, and ±g is the effective noise field arising from the 
coupling to the fluctuator. The average energy level sep- 
aration of the qubit is Bq in units where H = 1 and 
giiB = 1- Our notation is that of a spin qubit but the 
model obviously applies to any two-level quantum sys- 
tem. 

If g is in the z-direction, then we refer to dephasing 
noise, and if g is in the x-y plane, then we have relax- 
ational noise. 

We characterize the fluctuator as follows. Within a 
small time interval At the fluctuator state s' stays un- 
changed, but it changes to s in going to the next inter- 
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val with the conditional probabihty Wss' ■ The switch- 
ing probabiUties are parameterized as W |- = p + S 

and = p — S, and normalization requires that 

= 1— p— (5 and Ty__ = 1—p+S. To mimic telegraph 
noise we require p = jAt <C 1 and \S\ = \r]\At <C 1 with 
finite rates 7 and t]. For a fluctuator in thermal equilib- 
rium the stationary level population of the two fluctuator 
states is easily shown to be Ps = {p—sd)/2p= {-f—sri)/2^ 
by a standard detailed-balance argument [7]. 

A sequence of the noise simulated the model is shown 
in Fig.l. We emphasize that telegraph noise has a fi- 
nite correlation time, unlike white noise, for which exact 
methods are already available. Indeed, the telegraph se- 
quence is a Poisson process and the noise spectrum is 
Lorentzian. 

We wish to solve for the qubit density matrix p: given 
p{t = 0), find p at all later times. It has been shown that 
this decoherence problem is exactly solvable by meth- 
ods coming from the theory of stochastic differential 
equations. [8] The formal solution is given in terms of a 
complicated Laplace transform, which would be difficult 
to invert analytically and numerically unstable in most 
cases. In special cases the solution can be simplified and 
results have been presented when (1) 5a; = 5i/ 0: pure 
dephasing (T2) noise [4, 9-11] (and references therein) 
and (2) = 0: pure transverse or relaxational (Ti) noise 
[10, 11]. These special cases are important. However, 
there are clearly many situations in which both types 
of noise are present. For example, it has recently been 
shown that the character of the noise in flux qubits can be 
changed continuously from one type of noise to the other 
by changing the bias voltage [12]. The crossover region 
can only be described by the more general model, and 
the only available treatment is to use the phase-memory 
functional in the linear-coupling regime where g <C Bq 
[13]. In this paper, we show how to solve the general 
problem for a single qubit by a new algebraic; method. 
This greatly simplifies the special cases and makes possi- 
ble the presentation of results for all values of{gx,gy,gz}- 
The solution is exact and does not make perturbative ap- 
proximations. The main results are Eqs.(3)-(4) for the 
discrete-time formalism, Eqs.(6)-(7) for the continuous- 
time formalism, and the discussions thereafter. 

The remainder of this paper is organized as follows. We 
develop the exact solution of the problem by the general- 
ized transfer matrix method in section II, apply the the- 
ory to the case of bang-bang control in section III, and to 
the case of echo decay in section IV. We discuss in Sec.V 
possible generalization of the theory to the cases of com- 
bined white and telegraph noise, many fluctuators, and 
many qubits. Finally, Sec. VI is a summary of the work. 



II. TRANSFER MATRIX ENSEMBLE 
AVERAGED OVER TELEGRAPH NOISE 
SEQUENCES 

Our task is to solve for p{t = NAt) in the presence 
of a time-dependent Hamiltonian, averaged over all 2^ 
sequences with the appropriate probabilities. Formally 

p{t) = Us,Us,_, ■ ■ ■ u,,p{o)u!, ■ ■ ■ u!,_,u!,, 

where Us = exp {—iHgAt) is the evolution operator at 
noise level s for one time interval and the overbar indi- 
cates the averaging over the sequences. However, it is 
much more convenient to parameterize p by 

Pit) = 1^+ l^it) ■ ^- 

The Bloch vector n{t) takes on only real values and sat- 
isfies \n{t)\ < 1 (the equality holds for a pure state). / 
is the 2x2 unit matrix. Note that n{t) also gives the 
qubit magnetization as {a) = Tr(pa) = n{t). The time 
development is given by 

n (t) = Tn(0) = T,^T,^_, • • • T,, n(0) , (2) 

which defines the ensemble averaged transfer matrix T, 
and the matrix Tg is given by 

T, = exp [iAt {BoL^ + sg • L)] , 

where L^ y ^ arc the usual generators of S0{3): {Li)^^. = 
i^ijki where Cijk is the completely antisymmetric symbol. 
The Ts are easily calculated explicitly for arbitrary g. For 
noise that is uncorrelated between time intervals of equal 
length, the ensemble average in Eq.(2) can be performed 
within each interval, allowing exact solution of certain 
nois(^ models [14], For true telegraph noise the ensemble 
average has to be done in another way which we now 
describe. 

Let us define GJ| to be the 3x3 transfer matrix for 
an iV-step qubit evolution that starts at the fiuctuator 
state s' and end up at the state s, but ensemble averaged 
over all intermediate fiuctuator states. For = 1 no 
intermediate intervals are involved, so that Gf** = F** = 
Wsa'Tg'. Here Tg/ is the transfer matrix in the starting 
state s', and Wgs' signifies the conditional probability for 
the change to s immediately after the interval ends. By 
definition, we find that 

Gss' \ ^ T^ss" r^s' s' 
N — <-^JV-H 

s" 

which is already in the form of a matrix product, and 
by iteration we see that Gn = F^, where we defined an 
operator 

r = {1 - p - St3 + pTl - iST2 ) X 

exp{iAtBo-L + iAtg-LT3), (3) 
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where the PauU matrices {ri} act on the fluctuator state 
\s) (the eigenstate of T3). It is easily seen that F** = 
(s|r|s'). Finally the globally ensemble-averaged transfer 
matrix is given by 

T = 5^G?|W = (a^/|r^lv), (4) 

s,s' 

where \xf) = X^s=± is formally one of the eigen- 
state of Ti; and \if) = V^^g^^Ps\s) encodes the initial 
stationary level distribution of the fluctator. Note that 
the formal inner product is performed in the fluctuator 
level space, leaving a 3 x 3 matrix acting on the initial 
qubit vector n(0). 

In the hmit of Ai ^ 0, 

r ^ 1 + iAt{BoL, + g • Lra) 
-At'j + At{'yTi — ir]T2 — tjts) 
~exp(-AiP), (5) 

where we define 

P = 7 - iBoL^ - ig • Lra - 7T1 + i77T2 + 77x3. (6) 

Thus in the continuum time limit G = = exp{—tP) 
and so 

T = {xf\eM-tP)\if) (7) 

for t = TV At, and the problem reduces to the diagonaliza- 
tion of P which can be cast into a 6 x 6 matrix. Assuming 
that P is not defective (an assumption we have checked 
in the cases treated here), T can be decomposed as 

T = J2{xf\X){X\if)eM-^t), (8) 

X 

where [A) and (A[ arc the right and left cigen vectors of P 
with the eigenvalue A, normalized such that (A|A') = Sw. 
(Notice again the partial inner products with \if) and 
{xf\.) Each eigenvalue corresponds to a relaxation time 
of the system. Typically, the shorter times correspond 
to transients and the longest two times correspond to Ti 
and T2 - this can be verified by examining n(t)) in detail. 

Now let be the angle between g and the z— axis. For 
= 0, g = gz, Lz is conserved by P, so that we can 
use the quantum numbers m = 0, ±1 in place of L^, and 
diagonalize the 2x2 matrix in the fiuctuator spin space. 
The eigenvalues A are given by 

A = 7 — iBqui ± 1/7^ — g'^m? — 2igrim, m = 0, ±1. 

The eigenvectors, and finally the transfer matrix T can 
also be easily obtained. We shall not go into these details 
here. By inspection of Re(A), we see that the relaxation 
rate 1/Ti = in the nz{t) channel (to = 0), whereas 
the dephasing rate I/T2 = 7 — Rey/^'^ — g'^ ± 2igri in the 



Ux y channel. The result agrees with the results in Ref. 

[9]: 

The other case for which compact explicit expressions 
can be give is 9 = it/ 2 (so that g is perpendicular to Bq) 
and at the same time the switching rate imbalance 77 = 0. 
For simplicity let us assume that g is in x-direction and 
Bo is still in the z-dircction. Upon a rotation in the level 
space, we have UnU^ = T3 and UtzU^ = — n, where U 
is the SU(2) rotation about the y-axis by 90 degrees (in 
the level spin space). Under this transformation 

UPU"^ = 7 - 7T3 - iBqIz + igL^n. 

This matrix is the same as that of a Hamiltonian for the 
coupling of a spin-1 particle with angular momentum L 
and a spin- 1/2 particle with angular momentum S = r/2. 
Inspection reveals that the operator only mixes states 
whose z-components of L + S differ by 2. (This is seen 
by writing = {L+ + i-)/2 and = ti/2 = {S+ + 
S~)/2). As the result, the Hilbert space is divided into 
two invariant subspaces, spanned respectively by the two 
sets of states 

{|-1,-1/2),|1,-1/2),|0,1/2)}, 
{|-l,l/2),|0,-l/2),|l,l/2)}. 

Here in each basis state the first index refers to L^, and 
the second one to that of Sz = T3/2. The diagonal- 
ization can carried out separately in the two subspaces, 
and only 3x3 matrices are involved. The T-matrix 
is formally given by T = (z/|exp(-t[/PJ7t)|z/), where 
\zf) = (1,0)-^. Here we used the facts that = |a;/) 
for = and that U\xf) = \zf). Without going into 
details we mention that the eigenvalues of UPW (and 
thus of P) satisfy one of the two equations below: 

X^ + 2jX^ + {Bi + g^)X + 2B^oj = 0, (9) 
A^ + 47A2 + {Bl + g^+ Aj^)X + 2g^j = 0. (10) 

The result can also be obtained from the stochastic dif- 
ferential equation approach [10]. 

For general 6, the analysis (in terms of phase mem- 
ory function) in the literature is limited so far to the 
so-called linear-coupling regime g <^ Bq [13]. Our model 
is however exact for any coupling strength. All that we 
have to do is to diagonalize a 6 x 6 matrix to obtain the 
full qubit dynamics. In Fig. 2 we plot the evolution of 
the qubit vector n in the rotating frame defined by Bq. 
The parameters are Bq = 1, ^ = 0.1, g = 0.3, rj = 
and 9 = tt/A. We see from Fig. 2(a) that starting from 
n(0) = X, Ux decays with oscillations even in the rotating 
frame. On the other hand the ^-component is induced in 
the intermediate stage and actually decays more slowly 
than the x-component. We checked that this feature is 
visible in the strong coupling regime (7 < g), but is much 
weaker in the weak coupling regime (7 > g), and is com- 
pletely absent in the case ol "f g, the limit of white 



4 



•0.5 





n 

X 


\ (a) 


- - . n 




y 

n 

z 






y \/r;=^ ■ 




100 



20 



40 60 
t 



80 



100 



2: (Color online) The evolution of the qubit vector n = 
n„, Tiz) in the rotating frame, starting from (a) n(0) = £, 



FIG 

(jlx , lt,y , 

and (b) n(0) = z. The parameters are Bo = 1, 7 = 0.1, 
g — 0.3, »7 = 0, and 9 = 7r/4. The time t is in units of I/Bq. 
The same legend for the curves is used in (a) and (b). 




FIG. 3: (Color online) The relaxation rate l/Ti and dephas- 
ing rate I/T2 (in units of Bo) as functions of the angle 8 
between the static field Bo and the noise field g. I/T2 (trian- 
gle points) is I/T2 calculated from the perturbative formula 
in Eq. 11. The weak coupling case (a) has the parameters 
7 = 0.5 and g — 0.1. Here -q — (0.1) for solid (dashed) 
lines), while the strong coupling case (b) uses 7 = 0.1 and 
g — 0.3. Here rj — (0.05) for solid (dashed) lines. Rescal- 
ing of the rates is used to improve clarity. The perturbative 
formula breaks down completely for strong coupling. 



noise. This signifies the unique role of long time cor- 
relations. In Fig. 2(b) the qubit starts from n — z and 
decays slower than the x- and y-components induced in 
the intermediate stages, in agreement with the behavior 
in Fig.2(a). 

One should be aware that none of this qubit dynamics 
is available from Redfield theory, which applies only to 



the regime 7^ ^ 1. 

We can decide precisely the asymptotic decay rates in 
the different channels by matching nx^y^z{t) with the en- 
velope curves exp [—Re (A) i], where A are the numerical 
eigenvalues of P. The resulting decay rates l/Ti (in the 
n^-channel) and \/T2 (in the na;^j,-channel) as a func- 
tion of 9 are plotted in Fig. 3. First consider the case 
of r/ = (solid lines). Fig. 3(a) is in the weak-coupling 
regime, where we see a smooth change and a crossing of 
the two rates as 9 increases. This is similar to the case of 
uncorrelated noise[15]. The reason is that for 5 < 7, in a 
time scale determined by 1 / g many switches occur and as 
such the noise is essentially uncorrelated beyond a time 
scale of 1/g. Fig. 3(b) is in the strong-couphng regime, 
where we note the transition from the flat behavior of T2 
at small angles to a downturn at larger angles. This is 
because I/T2 is largely controlled by ijz, which decreases 
with increasing 6 and eventually falls into an "effective 
weak-coupling" regime gz < 1 for I/T2. We also observe 
that T2 = 2Ti at 9 = tt/2 for both g > 7 and g < 7, and 
is therefore a consequence of the model for all coupling 
strengths (provided that 77 = 0). This can be checked 
analytically from the eigenvalue equations in the partic- 
ular limits of 7 — !■ or f; — > 0. Second, we consider 
the effect of a nonzero rj (dashed lines). As compared 
to the 77 = case, we see that in the strong coupling 
case (b) the flat regime and the sharp transition of I/T2 
are smeared, and in both Fig. 3 (a) and Fig. 3(b) the decay 
rates become smaller. This is understandable from the 
fact that a nonzero rj amounts to a nonzero average of 
the noise field g and thus a decrease of the amount of 
fluctuating component. In particular, we have checked 
that 1/Ti = l/r2 = for ry = ±7, as one would have 
anticipated since the fluctuator stops at one of the two 
levels and does not switch at all. 

In the literature the decay rates at arbitrary 9 are avail- 
able only in the weak coupling cases. The perturbative 
results are often summarized as [161 



with 



i/r2 = 1/2T1 + 1/T4, 

l/T^ = cos^ 9S{0)/2, 
1/Ti = sin^ 0S{n)/2. 



(11) 



(12) 
(13) 



Here l/T^ is the dephasing rate as if in an effective 
z-direction random field with amplitude Bq cos 9 alone, 
S{uj) is the power spectrum of telegraph noise and ft = 
Bq. In view of the wide use of this formula, it is im- 
portant to check how well it holds for general coupling 
strengths. The comparison is made in Fig. 3 (for rj — 
only), where l/r2 calculated according to the above for- 
mula is denoted by l/Tj* (triangles) and thus can be 
compared to the actual I/T2. We see that in the weak 
(strong) coupling case of Fig. 3(a) [Fig. 3(b)] it agrees per- 
fectly with (deviates considerably from) our exact result 



5 



(solid lines). The formula cannot be used when the cou- 
pling of the source to the qubit is large compared to the 
inverse correlation time of the noise. This illustrates the 
importance of our exact results in the strong coupling 
cases. 

III. DYNAMICAL DECOUPLING BY 
BANG-BANG CONTROL 

The formalism can be easily adapted to echo decay 
measurements and the bang-bang control protocol where 
control pulses are applied. In such processes the transfer 
matrix can be formally written as 

T= (a;/|f exp[- / dt' P + i^, HU) ■ f) , 
Jo 

where T time-orders the operators, encoding the instan- 
taneous rotations with vectors caused by pulses 
at time {U}. These manipulations are important for com- 
parison with experiments such as that in Ref.[12]. 

Let us consider an open loop quantum control (or 
bang-bang control), in which a sequences of or tTj, 
pulses with fixed intervals are applied to reduce the de- 
coherence due to slowly low frequency noise. In the ideal 
case, we assume that each pulse is of zero width in time. 
The problem of ideal dynamical decoupling by bang-bang 
control has been exactly solved in the special cases of 
pure dephasing noise {9 = 0) and pure transverse noise 
(6* = 7r/2) [11]. With our method we can get the exact 
solution at arbitrary working points. 

In some cases ttx pulses would cause the anti-Zeno ef- 
fect in experiment [17], so that we consider iTy pulses. Af- 
ter a TTy pulse is applied, the Bloch vector is rotated by 
TT about the y-axis: n — > exp(i7rLj,)n. For a periodic 
sequence of tTj, pulses, the transfer matrix is explicitly 
given by 

T{Nt) = {xf\ [eM-Pr)eM^^Ly)f \if), (14) 

where r is the interval between two adjacent pulses, N is 
the number of pulses applied. By diagonalizing the oper- 
ator (again a 6 x 6 matrix) exp(— Pr) exp(i7rLy), we get 
the eigenvalues Xi=i^...fi and the candidate decay rates 
Ti — — ln]Ai]/r. In principle only one of the six F's 
control the long time asymptotic behavior in a specific 
channel. This rate can be decided precisely, as we used 
in Sec. II, by matching nx.y,z{t) with the envelope curve 
exp(-rjt). 

In Figs. 4 we plot the normalized decay rates Ti/Tid 
(in the n^-channel) and T2/T2d (in the n^^^j^-channel) as 
functions of the interval r between pulses. l/Ti and I/T2 
are the rates in the absence of the pulses. A normalized 
rate of 1 therefore corresponds to no suppression of de- 
coherence by the pulses. We set = 1, ?y = 0, 6* = 7r/4, 
and 7 = 0.1. In addition, g = 0.03 and g — i in Fig. 4(a) 
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FIG. 4; (a)The normalized rates for bang-bang control. 
The relaxation rate Ti/Tid (triangles) and the dephasing rate 
T2/r2d (squares) are plotted as functions of the interval r be- 
tween Hy pulses, with the parameters Bo = 1, 7 = 0.1, 
g = 0.03, r) = Q, and 9 — tt/4. (b) The same plot as (a) 
except that g — 3. 1/Ti and I/T2 are the rates in the absence 
of the pulses. Lines are drawn to guide the eye. 

and (b), respectively. We observe that as the interval r 
between iTy pulses decreases, the decay rates of both Uz- 
channel and n^j y-channel decrease. The reduction is sig- 
nificant as soon as t ~ 1/g. Moreover, in Fig. 4(b) where 

6 = tt/4: we observe oscillatory behavior of the relative 
rates similar to the case oid = tt/2 studied elsewhere. [11] 

IV. SPIN-ECHO DECAY 

In a spin-echo process, a pulse is applied half-way 
between two 7r^/2 pulses. The first 7r^/2 pulse rotate the 
initial quibt state to Uy, the second rotate Uy 

back to Uz right before the measurement. The tTx pulse 
reverses the ny-component, reducing the line-broadening 
associated with low frequency noise in the final measure- 
ment. The transfer matrix in the combined process is 
easily shown to be given by 

T{t) = {xf I exp(j7rL^/2) exp(-iP/2) exp(i7rL^) 

xexp{-tP/2)exp{iTTL^/2)\if), (15) 

where t is the time interval between the two nx/2 pulses. 

In Fig. 5 we plot the time dependence of the echo sig- 
nal for the parameters g — 0.8 (strong coupling case), 

7 = 0.1, and (a) 6* 0, (b) 61 = tt/A and (c) 9 = 7r/2. 
We can see that the echo signal shows steps. The 
steps in the special case of = have been observed 
experimentally, [18] as has been pointed out in Ref.[4]. 
Here we point out that steps should occur for all 9, and 
the period of the steps depends sensitively on 9. Indeed, 
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FIG. 5: The time dependence of echo signals for g = 0.8, 
7 = 0.1 and (a) 6 = 0(pure dephasing), (b) 6 = n/A and (c) 
6 = Tv/2. 
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FIG. 6: The time dependence of echo signal for g 
7 = 0.1. (a) gives the case 6 = 0, (b) is the case 9 



0.08 and 
7r/4 and 

(c) is the case 6 = vr/2, which corresponds to the experiment 
in Ref. [121. 



the longest period and the most pronounced steps occur 
for = tt/A. 

Figs. 6 are the same plots as Figs. 5 except that g — 0.08 
(weak coupling case). In this case no steps are ob- 
served. The signal decays exponentially. This is in agree- 
ment with the experiment at the optimal working point 
{6 = 7r/2)[12], but our results show that this is a general 
behavior for any values of 9 and the decay is faster for 
6 < tt/2 (which is reasonable since the z-component of 
the noise field increases as 9 decreases). The echo decay 
rate in the one fluctuator model is not reduced signifi- 
cantly as compared to free-induction decay, while in the 
experiment of Ref. [12, 18] the reduction is very signifi- 
cant. This indicates that a many- fluctuator model (plus 
possibly other types of noise) is needed. This will be 



discussed in the next section. 



V. GENERALIZATIONS 

The formalism developed so far can be extended to 
the case of combined noise sources: a white-noise field 
h and a two-level fluctuator. We assume that in a time 
sequence h does not change within an interval At, but 
is uncorrelated between different intervals. This enables 
us to do the ensemble average over h within each time 
interval, where we assume (hihj) - 
X, y, z. The effect of h enters F as 



(vil At)bi2, with i 



r 



(1 - p - (5r3 
exp[zAi(Bo - 



h) -L-hiAtg-Lra]. 



(16) 



Expanding the exponential function (up to the second 
order of h) and performing the ensemble average over 
h, and finally taking the limit At ^ (recalling that 
p — 7At and (5 = T^At) we again obtain F = cxp(— AiP) 
with a modified operator 

P = 7 + ^ v,L\l2 ~ zBoi^ 

-iS-'LT3-jTi+i7]T2 + r]T3. (17) 

Here we assumed that the contributions from higher or- 
der moments of h are of higher orders in At and can be 
ignored. As such the result is independent of the concrete 
form of the distribution function for h. 

The theory can also be generalized to many indepen- 
dent fluctuators described by the parameters {g„, jniVn} 
(n = 1, N). We assume that after each small time in- 
terval at most one of the fiuctuators could switch (joint 
switches occur at higher orders of At) . Since the fiuctua- 
tors are independent, it is straightforward to see that the 
desired P operator is now given by 

P = ^{ln - ign ■ Lt3„ - JnTln + irinT2rL + 77riT3„) 

n 

-iBoL,, (18) 

where Tin denotes the z-th PauH matrix operating on the 
71-th fiuctuator spin. This is a system with a spin-1 
coupled to many independent spin-l/2's. The problem 
is still solvable algebraically provided that [i2,P] — 0. 
In general, P is a (3 x 2^) x (3 x 2^) matrix. The fi- 
nal 3x3 T-matrix according to which the Bloch vector 
evolves is formally given by T = (a;/| exp(—tP)\if) with 
\xf) = Unlxfn) and \if) = n„|i/„), where \xfn) and |i/„) 
describe the states of the n-th fiuctuator. 

Finally the method can be generalized to the case of 
many qubits. Let D be the dimension of the Hilbert space 
of the qubits. Then there are — 1 Hermitian generators 
Xi of the transformation group SU{D) instead of the 
three Pauli matrices used in the single qubit case where 
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D ~ 2. The D X D density matrix may be expanded as 
a real linear combination of the Xi as in Eq.(2): 

D^-l 

p(t)=I/D+ J2 

The - l) X _ i) matrix T, that evolves the vec- 
tor a(t) forward in time when the joint state s of the A'' 
two-level fluctuators are given is now easily computed, 
and the steps in deriving the transfer matrix are formally 
identical to that in Sec. II. We see that the correspond- 
ing operator P wih be a {2^{D^ - 1)) x {2^{D^ - 1)) 
matrix, and its eigenvalues determine the (i?^ — l) re- 
laxation times in the system, one corresponding to each 
possible observable. The rapidly increasing dimension of 
the superoperator space will begin to give difficulties for 
numerical calculations even at moderate values of A'' and 
D. 

VI. SUMMARY 

We have developed an exact transfer matrix method 
to solve the problem of qubit decoherence caused by one 
fluctuator, and have applied the theory to free qubit de- 
cay, bang-bang control and spin-echo decay. We have 
reproduced all known exact results and have shown that 
the perturbative limits are correct. The method is rela- 
tively straightforward to apply, since it is completely al- 
gebraic, and the qubit decoherence rates are determined 
by the eigenvalues of a transfer matrix. 

The formalism is also more powerful than previous ex- 
act methods in that it can be applied to a qubit subject to 
an arbitrary combination of dephasing and relaxational 
telegraph noise. The typical way of combinnig the two 
types of noise is according to the perturbative formula 
I/T2 = l/2Ti + l/T^. We have shown that this formula 
breaks down when the coupling is strong. The algebraic 
method can get qubit dynamics on time scales short com- 
pared with the environmental correlation times, which is 
often of experimental interest. 

The method can be generalized to the case of many 
noise sources and multiple qubits, though the size of the 
matrices grows rapidly. 
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